% clear all

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Date: March 2, 2018
%
% Program: Calls the elasticity decomposition function and saves the
%          results. This function decomposes aggregate elasticities given
%          shocks at the firm and sector levels. 
%
% NEW: preapred to run with CES coding, so also change notation. Further,
% call up data in master file rather than in this step separately.
%
% Updated: 8/2/2019 by IM/JDG
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
global alphaT rhoT lambdaT etaT shockT rho sigma sigmaT lambda eta ...
    tol maxit D_hat_n Xi_hat_mnjf Xi_hat_mnj a_hat_f a_hat ...
    varphi sL_n sPi_n sD_n N J FI_J france oneminus_alpha_vector ...
    labour_share_PWT countrysecD firmsecD_sorted start start_sorted ...
    sum_by_country_dummy ISO alpha_WIOD_j alpha_WIOD_francej sector_algorithm 


% Compute weights
% Firm-level VA shares
% eval([strcat('load MAT/rho',num2str(rhoT),'eta',num2str(etaT),'lambda',num2str(lambdaT),...
%     '/Step_4_setup_alpha',num2str(alphaT),...
%      '_rho',rhoT,'.mat')]);
% X_mnjf_0 = firmsecD_sorted*X_mnj(start(france):start(france+1)-1,:);
% rho_f=firmsecD_sorted*rho;
% pi_l_f_repmat0=repmat(pi_l_f, [1 N]); 
% rho_fn=repmat(rho_f,[1 N]);
% VA_fnj_0 =((1+pi_l_f_repmat0.*(rho_fn-1))./rho_fn).*X_mnjf_0;
% VA_fj_0 = nansum(VA_fnj_0,2);
% 
% % Sector-level VA shares
% rho_mnj_repmat=repmat(rho,[N N]);
% pi_l_mnj_0 = repmat(pi_l,[1 N]);
% VA_mnj_0=((1+pi_l_mnj_0.*(rho_mnj_repmat-1))./rho_mnj_repmat).*X_mnj;
% VA_mj_0 = nansum(VA_mnj_0,2);

%N = size(VA_hat_n,1);
%J = size(TFP_hat_cpi_france_j,1);

% French Firm results

[ e_Y_FRf_cpi, m_e_f_cpi, cov_e_f_cpi, check_f_cpi,...
     e_Y_FRf_doubledef,m_e_f_doubledef,cov_e_f_doubledef,check_f_doubledef,...
     e_Y_FRf_def, m_e_f_def, cov_e_f_def, check_f_def]...
     = elasticity_decomp_firmCES_fun(delta,france,RGDP_cpi_hat_fj,RGDP_cpi_hat_n,...
        VA_fj_0,RGDP_def_hat_fj,RGDP_def_hat_n,RGDP_doubledef_hat_fj,RGDP_doubledef_hat_n);

% Sector France Results

s1 = start(france);
s2 = start(france+1)-1;
[e_Y_FRs_cpi,m_e_mjFRs_cpi,cov_e_mjFRs_cpi,check_mjFRs_cpi,...
    e_Y_FRs_doubledef,m_e_mjFRs_doubledef,cov_e_mjFRs_doubledef,check_mjFRs_doubledef,...
    e_Y_FRs_def,m_e_mjFRs_def,cov_e_mjFRs_def,check_mjFRs_def]...
    = elasticity_decomp_sectorFRCES_fun(delta,france,s1,s2,VA_hat_mj,P_hat_n,...
    RGDP_cpi_hat_n,VA_mj_0,GDP_deflator_n,RGDP_def_hat_n,DoubleDeflation_deflator_n,RGDP_doubledef_hat_n);
 
% Sector results

[e_Y_s_cpi,m_e_mj_cpi,cov_e_mj_cpi,check_mj_cpi,...
    e_Y_s_doubledef,m_e_mj_doubledef,cov_e_mj_doubledef,check_mj_doubledef,...
    e_Y_s_def,m_e_mj_def,cov_e_mj_def,check_mj_def] ...
    = elasticity_decomp_sectorCES_fun(delta,VA_hat_mj,...
    P_hat_n,RGDP_cpi_hat_n,VA_mj_0,GDP_deflator_n,RGDP_def_hat_n,DoubleDeflation_deflator_n,RGDP_doubledef_hat_n);

display('Computation of elasticities completed');